Aula09 - Introdução à Ciência dos Dados

Armazenamento de números em computadores

Qual seria o resultado da seguinte subtração no sistema decimal?

(0.3 - 0.2) - 0.1

## [1] -2.775558e-17

Alguma ideia sobre a causa do resultado?

Ponto Flutuante Binário

A representação em ponto flutante é uma variação binária (base-dois) da notação científica, sendo esta a representação utilizada por computadores para armazenar números.

Por exemplo, podemos escrever o número \(0,0006926\) com quatro dígitos significativos em notação científica como \(6,926 × 10^{-4}\). Esta representação pode representar qualquer valor entre \(0,00069255\) e \(0,00069265\).

As representações em ponto flutuante em computadores são semelhantes, exceto que uma potência de 2 é usada em vez de uma potência de 10, e a fração é escrita em notação binária.

O número \(0,0006926\) é escrito como \(1,0112 \times 2^{-11}\) se a precisão de quatro dígitos binários for usada. \(1,0112\) representa a mantissa e \(11\) o expoente.

No entanto, \(6,926 × 10^{−4}\) e \(1,0112 × 2^{−11}\) não são idênticos.

Quatro dígitos binários fornecem menos precisão que quatro dígitos decimais. Uma faixa de valores entre 0,000641 a 0,000702 obteria a mesma representação com precisão de quatro dígitos binários.

De fato, \(6,926 × 10^{−4}\) não pode ser representado exatamente em notação binária em um número finito de dígitos.

O problema é semelhante a tentar representar \(1/3\) como um decimal: \(0,3333\) é apenas uma aproximação.

.Machine é uma variável que contém informações sobre as características numéricas da máquina em que a linguagem R está sendo executada, como o maior double ou inteiro e a precisão da máquina.

.Machine

## $double.eps
## [1] 2.220446e-16
## 
## $double.neg.eps
## [1] 1.110223e-16
## 
## $double.xmin
## [1] 2.225074e-308
## 
## $double.xmax
## [1] 1.797693e+308
## 
## $double.base
## [1] 2
## 
## $double.digits
## [1] 53
## 
## $double.rounding
## [1] 5
## 
## $double.guard
## [1] 0
## 
## $double.ulp.digits
## [1] -52
## 
## $double.neg.ulp.digits
## [1] -53
## 
## $double.exponent
## [1] 11
## 
## $double.min.exp
## [1] -1022
## 
## $double.max.exp
## [1] 1024
## 
## $integer.max
## [1] 2147483647
## 
## $sizeof.long
## [1] 8
## 
## $sizeof.longlong
## [1] 8
## 
## $sizeof.longdouble
## [1] 8
## 
## $sizeof.pointer
## [1] 8
## 
## $sizeof.time_t
## [1] 8

• A constante double.eps representa o menor número positivo que pode ser representado no tipo double antes que a precisão seja perdida, é também conhecida como machine epsilon

abs((0.3 - 0.2) - 0.1) < .Machine$double.eps

## [1] TRUE

• Podemos ver a precisão dupla (daí o uso de double) da linguagem R (e náo apenas da linguagem R) é de 53 bits (dígitos binários), o que equivale a de 15 ou 16 dígitos decimais.

• Assim, a aritmética com números inteiros será exata para valores entre \(−(2^{53} − 1)\) e \(2^{53} − 1\) (aproximadamente \(−10^{16}\) a \(10^{16}\)), mas assim que começamos a usar números fora desse intervalo, ou frações, podemos perder precisão devido a erros de arredondamento.

• Por exemplo, \(1.1\) não possui uma expansão binária finita, portanto, em precisão dupla (double precision), sua expansão binária é arredondada para \(1,00011001100\ldots001\) , com um erro de aproximadamente \(2^{-53}\) .

• O erro de arredondamento tende a se acumular na maioria dos cálculos, então geralmente uma longa série de cálculos resultará em erros maiores do que uma curta.

• Algumas operações são particularmente propensas a erros de arredondamento: por exemplo, subtração de dois números quase iguais, ou (equivalentemente) adição de dois números com quase a mesma magnitude, mas sinais opostos.

• Isso se ocorre por que os bits principais nas expansões binárias de números quase iguais serão correspondentes, assim serão cancelados na subtração, e o resultado vai depender do que está armazenado nos bits posteriores.

Exemplo:

Considere a seguinte fórmula padrão do estimador da variância amostral de uma amostra \(x_1,\ldots,x_n\):

\[ s^2 = \frac{1}{n-1} \sum_{i=1}^n (x_i - \bar{x})^2 \]

sendo \(\bar{x} = 1/n\sum x_i\). Como sabemos, em R, \(s^2\) é obtido por var() e \(\bar{x}\) por mean(). Por exemplo

x <- 1:11
mean(x)

## [1] 6

var(x)

## [1] 11

sum((x - mean(x))^2)/10

## [1] 11

Como este estimador requer o cálculo de \(\bar{x}\) primeiro e a soma de desvios quadráticos em seguinda, isto requer que todos os valores \(x_i\) sejam mantidos na memória.

Não muito tempo atrás, a memória era tão cara que era vantajoso reescrever a fórmula como:

\[ s^2 = \frac{1}{n-1} \sum_{i=1}^n (x_i^2 - n\bar{x}^2) \]

Este estimador é conhecido como “one-pass formula” (fórmula de passagem única) por que avaliamos cada valor \(x_i\) apenas uma vez, e acumulamos as somas de \(x_i\) e de \(x_i^2\). Ele fornece a resposta correta em nosso exemplo:

(sum(x^2) - 11 * mean(x)^2) / 10

## [1] 11

Entretanto, vejo o que ocorre se um grande valor A a cada \(x_i\). A soma \(\sum_{i=1}^n x_i^2\) aumenta em aproximadamente \(nA^2\), assim como \(n\bar{x}^2\), ou seja, fazemos a subtração de números muito próximos.

Isso não altera a variância, mas fornece as condições para uma “perda catastrófica de precisão” quando calculamos a diferença:

A <- 1.e10
x <- 1:11 + A
var(x)

## [1] 11

(sum(x^2) - 11*mean(x)^2) / 10

## [1] -13107.2

Tipos de Dados Tradicionais em Finanças

Dados em seção cruzada ou transversal (cross-section)

• Dados em secção cruzada ou transversal (cross-section):

• as observações coletadas em um determinado ponto no tempo sobre diferentes unidades de observação. Cada observação representa uma única unidade de observação, como indivíduos, empresas ou países, e fornece informações sobre várias variáveis dessa unidade em um determinado momento. A ordem das observações não é relevante.

• Um conjunto de dados coletado pela amostragem de uma população em um determinado ponto no tempo.

• Exemplo: dados sobre diversas variáveis de 200 empresas para o ano de

Exemplo:

O conjunto de dados Hmda é um conjunto de dados coletados em Boston/EUA, especificamente, é uma seção cruzada entre 1997-1998.

O Home Mortgage Disclosure Act foi promulgado para monitorar o acesso de minorias de baixa renda ao mercado de hipotecas. Os dados recolhidos para esse objetivo mostraram que as minorias têm duas vezes mais chances de ter uma hipoteca negada do que os brancos. No entanto, as variáveis correlacionadas com raça e capacidade de crédito foram omitidas desses dados, tornando impossível qualquer conclusão sobre o papel da raça nos empréstimos hipotecários. O Federal Reserve Bank de Boston coletou variáveis adicionais importantes para a decisão de empréstimo hipotecário e descobriu que a raça continuou a desempenhar um papel importante, embora significativamente diminuído, na decisão de conceder um hipoteca.

Dicionário dos Dados

• número de observações: 2381

• observação: indiíviduos

• fonte: Federal Reserve Bank of Boston

dir: razão entre o valor das prestações e a renda total

hir: razão entre despesas com moradia e renda

lvr: razão entre o tamanho do empréstimo e o valor estimado da propriedade

ccs: pontuação de crédito ao consumidor de 1 a 6 (um valor baixo sendo uma boa pontuação)

mcs: pontuação de crédito hipotecário de 1 a 4 (um valor baixo sendo uma boa pontuação)

pbcr: histórico de crédito público ruim?

dmi: seguro hipotecário negado?

self: autônomo?

single: o requerente é solteiro?

uria: taxa de desemprego de Massachusetts em 1989 na indústria do requerente

condominium: a unidade é um condomínio?

black: o requerente é negro?

deny: pedido de hipoteca negado?

data(Hdma, package = "Ecdat")
head(Hdma)

##     dir   hir       lvr ccs mcs pbcr dmi self single uria comdominiom black
## 1 0.221 0.221 0.8000000   5   2   no  no   no     no  3.9           0    no
## 2 0.265 0.265 0.9218750   2   2   no  no   no    yes  3.2           0    no
## 3 0.372 0.248 0.9203980   1   2   no  no   no     no  3.2           0    no
## 4 0.320 0.250 0.8604651   1   2   no  no   no     no  4.3           0    no
## 5 0.360 0.350 0.6000000   1   1   no  no   no     no  3.2           0    no
## 6 0.240 0.170 0.5105263   1   1   no  no   no     no  3.9           0    no
##   deny
## 1   no
## 2   no
## 3   no
## 4   no
## 5   no
## 6   no

Uma operação comum qunado a variável respota é binária, é transformar sim e não para 1 e o, respectivamente:

hdma <- Hdma |> 
          mutate(deny_int = as.integer(ifelse(deny == "yes", 1, 0))) |>
          glimpse()

## Rows: 2,381
## Columns: 14
## $ dir         <dbl> 0.221, 0.265, 0.372, 0.320, 0.360, 0.240, 0.350, 0.280, 0.…
## $ hir         <dbl> 0.221, 0.265, 0.248, 0.250, 0.350, 0.170, 0.290, 0.220, 0.…
## $ lvr         <dbl> 0.8000000, 0.9218750, 0.9203980, 0.8604651, 0.6000000, 0.5…
## $ ccs         <dbl> 5, 2, 1, 1, 1, 1, 1, 2, 2, 2, 1, 1, 1, 1, 1, 2, 1, 2, 2, 2…
## $ mcs         <dbl> 2, 2, 2, 2, 1, 1, 2, 2, 2, 1, 2, 2, 2, 1, 1, 1, 1, 1, 2, 2…
## $ pbcr        <fct> no, no, no, no, no, no, no, no, no, no, no, no, no, no, no…
## $ dmi         <fct> no, no, no, no, no, no, no, no, yes, no, no, no, yes, no, …
## $ self        <fct> no, no, no, no, no, no, no, no, no, no, no, no, no, no, no…
## $ single      <fct> no, yes, no, no, no, no, yes, no, no, yes, yes, yes, no, n…
## $ uria        <dbl> 3.9, 3.2, 3.2, 4.3, 3.2, 3.9, 3.9, 1.8, 3.1, 3.9, 3.1, 4.3…
## $ comdominiom <dbl> 0, 0, 0, 0, 0, 0, 1, 0, 0, 0, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0…
## $ black       <fct> no, no, no, no, no, no, no, no, no, no, no, no, no, no, no…
## $ deny        <fct> no, no, no, no, no, no, no, no, yes, no, no, no, yes, no, …
## $ deny_int    <int> 0, 0, 0, 0, 0, 0, 0, 0, 1, 0, 0, 0, 1, 0, 0, 0, 0, 0, 0, 0…

Análise Exploratória dos Dados

summary(hdma)

##       dir              hir              lvr              ccs       
##  Min.   :0.0000   Min.   :0.0000   Min.   :0.0200   Min.   :1.000  
##  1st Qu.:0.2800   1st Qu.:0.2140   1st Qu.:0.6528   1st Qu.:1.000  
##  Median :0.3300   Median :0.2600   Median :0.7794   Median :1.000  
##  Mean   :0.3308   Mean   :0.2553   Mean   :0.7378   Mean   :2.116  
##  3rd Qu.:0.3700   3rd Qu.:0.2988   3rd Qu.:0.8684   3rd Qu.:2.000  
##  Max.   :3.0000   Max.   :3.0000   Max.   :1.9500   Max.   :6.000  
##       mcs          pbcr       dmi         self      single          uria       
##  Min.   :1.000   no  :2205   no :2333   no  :2103   no :1444   Min.   : 1.800  
##  1st Qu.:1.000   yes : 175   yes:  48   yes : 277   yes: 937   1st Qu.: 3.100  
##  Median :2.000   NA's:   1              NA's:   1              Median : 3.200  
##  Mean   :1.721                                                 Mean   : 3.774  
##  3rd Qu.:2.000                                                 3rd Qu.: 3.900  
##  Max.   :4.000                                                 Max.   :10.600  
##   comdominiom     black       deny         deny_int     
##  Min.   :0.0000   no :2042   no :2096   Min.   :0.0000  
##  1st Qu.:0.0000   yes: 339   yes: 285   1st Qu.:0.0000  
##  Median :0.0000                         Median :0.0000  
##  Mean   :0.2881                         Mean   :0.1197  
##  3rd Qu.:1.0000                         3rd Qu.:0.0000  
##  Max.   :1.0000                         Max.   :1.0000

Usando um spineplot:

spineplot(deny ~ dir, data = hdma)

Spine plots podem ser vistos como uma generalização de gráficos de barras empilhadas e também de histogramas.

Podemo verificar que a obtenção de uma negativa aumenta à medida que a razão prestaçaõ/renda aumenta e que a relação não parece ser linear.

Obs. Note que a variável resposta (à esquerda de ~) precisa ser um fator e que a função categoriza a variável numérica dir que representa a razão prestaçaõ/renda.

Dados em seções cruzadas combinadas (pooled cross-section)

• Uma configuração de dados em que seções cruzadas independentes,
  geralmente coletadas em diferentes pontos do tempo, são combinadas para produzir um único conjunto de dados.

• Exemplo: dados sobre diversas variáveis de 200 empresas diferentes para os anos de 2009 e 2010.

Exemplo

Uma seção cruzada utilizada para precificar casas com 321 observações sobre 19 variáveis

year: 1978, 1981

age: idadde da casa

nbh: vizinhança, 1-6

cbd: dist. para cento. ônibus. dstrct, ft.

inst: dist. para interestadual, ft.

price: preço de venda

rooms: número de quartos na casa

area: metragem quadrada da casa

land: metragem quadrada do terreno

baths: número de banheiros

dist: dist. de casa para incin., ft

y81: = 1 se year = 1981

data(hprice3, package = "wooldridge")

hprice3_slice <- dplyr::slice(hprice3, 177:182) |> 
                  dplyr::select(year, age, agesq, nbh, cbd, inst, linst, price, 
                                rooms, area, land, baths, dist) |> 
                            data.frame()

head(hprice3_slice)

##   year age agesq nbh   cbd  inst  linst  price rooms area  land baths  dist
## 1 1978   0     0   0 14000 13000 9.4727 119900     7 2352 31331     3 20500
## 2 1978   0     0   0 14000 13000 9.4727 142500     7 2492 27248     3 20300
## 3 1978   0     0   0 14000 13000 9.4727  84211     7 2756 27205     3 20200
## 4 1981  81  6561   4  4000  1000 6.9078  49000     6 1554  6790     1 11800
## 5 1981  71  5041   4  3000  2000 7.6009  52000     5 1575  3485     1 10100
## 6 1981  31   961   4  3000  2000 7.6009  68000     6 3304 18731     2 10200

Análise Exploratória dos Dados

Utilizando o pacote correlationfunnel para identificar as variáveis mais correlacionadas com price:

hprice3_slice |>
    glimpse() |>
    select(-year) |>
    correlate(target = price) |>
    plot_correlation_funnel()

## Rows: 6
## Columns: 13
## $ year  <int> 1978, 1978, 1978, 1981, 1981, 1981
## $ age   <int> 0, 0, 0, 81, 71, 31
## $ agesq <dbl> 0, 0, 0, 6561, 5041, 961
## $ nbh   <int> 0, 0, 0, 4, 4, 4
## $ cbd   <dbl> 14000, 14000, 14000, 4000, 3000, 3000
## $ inst  <dbl> 13000, 13000, 13000, 1000, 2000, 2000
## $ linst <dbl> 9.4727, 9.4727, 9.4727, 6.9078, 7.6009, 7.6009
## $ price <dbl> 119900, 142500, 84211, 49000, 52000, 68000
## $ rooms <int> 7, 7, 7, 6, 5, 6
## $ area  <int> 2352, 2492, 2756, 1554, 1575, 3304
## $ land  <dbl> 31331, 27248, 27205, 6790, 3485, 18731
## $ baths <int> 3, 3, 3, 1, 1, 2
## $ dist  <dbl> 20500, 20300, 20200, 11800, 10100, 10200

Utilizando o pacote GGally para a mesma finalidade:

ggpairs(hprice3_slice, columns = c("price", "inst", "baths", "nbh"))

Dados em séries temporais

• são dados coletadas em ao longo do tempo (em geral, em intervalos
  regulares) sobre uma ou mais variáveis de uma unidade.

• Exemplo: Lucro líquido, dívida, etc para uma empresa ao longo de 25 anos (ou trimestres, meses, dias…).

Exemplo

Dados diários sobre o preço de fechamento das 30 ações que compõe o índice Dow Jones.

data(DowJones30, package = "fBasics")
head(DowJones30, 5)

##    X.Y..m..d   AA  AXP     T    BA  CAT    C   KO    DD    EK   XOM   GE    GM
## 1 1990-12-31 5.92 4.70 14.67 21.77 9.30 1.87 9.88 13.15 23.02 12.09 4.63 19.72
## 2 1991-01-02 5.92 4.70 14.67 21.77 9.30 1.87 9.88 13.15 23.02 12.09 4.63 19.72
## 3 1991-01-03 5.88 4.73 14.73 21.71 9.15 1.87 9.60 12.83 22.95 12.15 4.53 19.57
## 4 1991-01-04 5.76 4.73 14.79 22.50 9.00 1.92 9.85 13.06 22.74 12.27 4.47 19.00
## 5 1991-01-07 5.72 4.58 14.67 21.65 8.87 1.89 9.69 12.83 22.11 12.12 4.40 18.35
##    HWP   HD  HON INTC   IBM    IP  JPM  JNJ  MCD   MRK MSFT   MMM    MO    PG
## 1 3.91 2.80 5.70 1.20 27.86 20.45 2.67 7.31 7.08 11.83 2.08 31.41 10.55 20.91
## 2 3.91 2.80 5.70 1.20 27.86 20.45 2.67 7.31 7.08 11.83 2.08 31.41 10.55 20.91
## 3 3.91 2.74 5.67 1.20 27.95 20.21 2.67 7.23 6.99 11.45 2.09 30.90 10.27 20.70
## 4 3.85 2.75 5.65 1.20 27.86 20.50 2.69 7.18 7.02 11.27 2.11 30.80 10.24 20.48
## 5 3.88 2.65 5.50 1.19 27.39 20.45 2.58 6.96 6.80 10.97 2.08 30.43 10.08 20.05
##    SBC  UTX  WMT  DIS
## 1 9.75 9.31 7.05 7.86
## 2 9.75 9.31 7.05 7.86
## 3 9.64 9.19 7.05 7.82
## 4 9.56 9.24 6.97 7.79
## 5 9.40 8.95 6.82 7.58

Importando e Estruturando Dados de Séries Temporais

Caso 1: Arquivo de Dados náo Contém Datas

O arquivo unemployment.csv contém 63 dados mensais sobre a taxa de desemprego civil dos EUA como porcentagem da força de trabalho para os anos 2012-2017, obtidos da base da dados do Federal Reserve Bank of Saint Louis.

1. Importe os dados utlizando a função read_csv do pacote readr e nomeie o objeto como dados.

path1 <- "dados/unemployment.csv"
dados <- readr::read_csv(path1)

2. Transforme dados para a classe ts indicando que os dados são mensais e o primeiro período é janeiro de 2012 e nomeie o objeto como desemprego_ts.

desemprego_ts <- ts(dados, start = c(2012, 1), frequency = 12)
desemprego_ts

##      Jan Feb Mar Apr May Jun Jul Aug Sep Oct Nov Dec
## 2012 8.8 8.7 8.4 7.7 7.9 8.4 8.6 8.2 7.6 7.5 7.4 7.6
## 2013 8.5 8.1 7.6 7.1 7.3 7.8 7.7 7.3 7.0 7.0 6.6 6.5
## 2014 7.0 7.0 6.8 5.9 6.1 6.3 6.5 6.3 5.7 5.5 5.5 5.4
## 2015 6.1 5.8 5.6 5.1 5.3 5.5 5.6 5.2 4.9 4.8 4.8 4.8
## 2016 5.3 5.2 5.1 4.7 4.5 5.1 5.1 5.0 4.8 4.7 4.4 4.5
## 2017 5.1 4.9 4.6

3. Faça um gráfico do objeto desemprego_ts

plot(desemprego_ts, col = "blue")
grid()

3. Utilizando a função as.xts do pacote xts, converta o objeto desemprego_ts para um objeto da classe xts e nomeie este objeto como desemprego_xts. Faça um gráfico do objeto desemprego_xts com plot(desemprego_xts).

desemprego_xts <- as.xts(desemprego_ts)
plot(desemprego_xts)

5. Observando o gráfico da série temporal, descreva se existem padrões relevantes.

A inspeção visual do gráfico da série indica uma clara tendência de queda da taxa de desemprego. Além disso, vê-se um padrão sazonal mas irregular.

Caso 2: Arquivo de Dados Contém Datas

O arquivo ipeadata[06-04-2021-04-16].csv, obtido no site ipeadata contém dados diários sobre a taxa de câmbio real/dólar, sendo que a primeira coluna contém os dias.

Podemos definir que a primeira coluna são dias usando o argumento col_type() da funcào read_csv2 do pacote readr

1. Importe o arquivo ipeadata[06-04-2021-04-16].csv e nomeie o objeto importado como cambio_csv

file_csv <- 'dados/ipeadata[06-04-2021-04-16].csv'

cambio_csv <- read_csv2(file_csv,
                    na = c(" "),
                    skip = 1,
                    locale = locale(decimal_mark = ","),
                    col_names = c("dia", "real_dolar"),
                    col_types = cols(
                    dia = col_date(format = "%d/%m/%Y"),
                    real_dolar = col_double())
                    )

dplyr::glimpse(cambio_csv)

## Rows: 13,243
## Columns: 3
## $ dia        <date> 1985-01-02, 1985-01-03, 1985-01-04, 1985-01-05, 1985-01-06…
## $ real_dolar <dbl> 1.152000e-09, 1.152000e-09, 1.152000e-09, NA, NA, 1.173818e…
## $ X3         <lgl> NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA,…

2. Elimine a coluna X3 e retire os dados faltantes (NA):

cambio_csv <- cambio_csv %>% select(dia, real_dolar) %>% drop_na()
summary(cambio_csv)

##       dia               real_dolar    
##  Min.   :1985-01-02   Min.   :0.0000  
##  1st Qu.:1994-02-07   1st Qu.:0.1858  
##  Median :2003-02-26   Median :1.8190  
##  Mean   :2003-02-24   Mean   :1.7745  
##  3rd Qu.:2012-03-19   3rd Qu.:2.6707  
##  Max.   :2021-04-05   Max.   :5.9366

3. Gráfico de linha usando o ggplot2

cambio_csv %>%
  ggplot(aes(x = dia, y = real_dolar)) +
  geom_line(color = "#69b3a2") +
  labs(
    title = "Taxa de Câmbio diária R$/US$",
    subtitle = "comercial - compra - média",
    caption = "Fonte: ipeadata - Bacen/SGS",
    x = "dia",
    y = "Taxa de Câmbio (R$/US$)"
  ) +
  theme_minimal()

4. Gráfico interativo usando o pacote plotly

interativo <- cambio_csv %>%
  ggplot(aes(x = dia, y = real_dolar)) +
  geom_area(fill = "#69b3a2", alpha = 0.5) +
  geom_line(color = "#69b3a2") +
  xlab("dia") +
  ylab("Taxa de Câmbio Diária (R$/US$)") +
  theme_minimal()

ggplotly(interativo)

Dados em painel (longitudinais)

• Estrutura de dados construída a partir de seções cruzadas repetidas ao longo do tempo.

• Com um painel balanceado, as mesmas unidades de observação aparecem em cada período de tempo. Com um painel desbalanceado, algumas unidades não aparecem em cada período de tempo, muitas vezes devido ao atrito.

• Exemplo: dados sobre diversas variáveis de 50 empresas durante 5 anos
  (“micropainel”), ou sobre 12 países durante 25 anos (“macropainel”).

Exemplo 1

Os dados de Grunfeld (1958 ) compreendeem observações referentes a \(11\) grandes empresas estadunidenses para os anos de 1935-1954 em relação a 4 variáveis:

• investimento bruto real (invest),
• valor real da empresa (value) e,
• valor real do estoque de capital (capital)
• ano (year)

Originalmente empregados em um estudo sobre os determinantes do investimento corporativo em uma tese de doutorado da Universidade de Chicago, esses dados têm sido um clássico dos livros didáticos desde a década de 1970.

data(Grunfeld, package = "AER")
head(Grunfeld, 60)

##    invest  value capital             firm year
## 1   317.6 3078.5     2.8   General Motors 1935
## 2   391.8 4661.7    52.6   General Motors 1936
## 3   410.6 5387.1   156.9   General Motors 1937
## 4   257.7 2792.2   209.2   General Motors 1938
## 5   330.8 4313.2   203.4   General Motors 1939
## 6   461.2 4643.9   207.2   General Motors 1940
## 7   512.0 4551.2   255.2   General Motors 1941
## 8   448.0 3244.1   303.7   General Motors 1942
## 9   499.6 4053.7   264.1   General Motors 1943
## 10  547.5 4379.3   201.6   General Motors 1944
## 11  561.2 4840.9   265.0   General Motors 1945
## 12  688.1 4900.9   402.2   General Motors 1946
## 13  568.9 3526.5   761.5   General Motors 1947
## 14  529.2 3254.7   922.4   General Motors 1948
## 15  555.1 3700.2  1020.1   General Motors 1949
## 16  642.9 3755.6  1099.0   General Motors 1950
## 17  755.9 4833.0  1207.7   General Motors 1951
## 18  891.2 4924.9  1430.5   General Motors 1952
## 19 1304.4 6241.7  1777.3   General Motors 1953
## 20 1486.7 5593.6  2226.3   General Motors 1954
## 21  209.9 1362.4    53.8         US Steel 1935
## 22  355.3 1807.1    50.5         US Steel 1936
## 23  469.9 2676.3   118.1         US Steel 1937
## 24  262.3 1801.9   260.2         US Steel 1938
## 25  230.4 1957.3   312.7         US Steel 1939
## 26  361.6 2202.9   254.2         US Steel 1940
## 27  472.8 2380.5   261.4         US Steel 1941
## 28  445.6 2168.6   298.7         US Steel 1942
## 29  361.6 1985.1   301.8         US Steel 1943
## 30  288.2 1813.9   279.1         US Steel 1944
## 31  258.7 1850.2   213.8         US Steel 1945
## 32  420.3 2067.7   132.6         US Steel 1946
## 33  420.5 1796.7   264.8         US Steel 1947
## 34  494.5 1625.8   306.9         US Steel 1948
## 35  405.1 1667.0   351.1         US Steel 1949
## 36  418.8 1677.4   357.8         US Steel 1950
## 37  588.2 2289.5   342.1         US Steel 1951
## 38  645.5 2159.4   444.2         US Steel 1952
## 39  641.0 2031.3   623.6         US Steel 1953
## 40  459.3 2115.5   669.7         US Steel 1954
## 41   33.1 1170.6    97.8 General Electric 1935
## 42   45.0 2015.8   104.4 General Electric 1936
## 43   77.2 2803.3   118.0 General Electric 1937
## 44   44.6 2039.7   156.2 General Electric 1938
## 45   48.1 2256.2   172.6 General Electric 1939
## 46   74.4 2132.2   186.6 General Electric 1940
## 47  113.0 1834.1   220.9 General Electric 1941
## 48   91.9 1588.0   287.8 General Electric 1942
## 49   61.3 1749.4   319.9 General Electric 1943
## 50   56.8 1687.2   321.3 General Electric 1944
## 51   93.6 2007.7   319.6 General Electric 1945
## 52  159.9 2208.3   346.0 General Electric 1946
## 53  147.2 1656.7   456.4 General Electric 1947
## 54  146.3 1604.4   543.4 General Electric 1948
## 55   98.3 1431.8   618.3 General Electric 1949
## 56   93.5 1610.5   647.4 General Electric 1950
## 57  135.2 1819.4   671.3 General Electric 1951
## 58  157.3 2079.7   726.1 General Electric 1952
## 59  179.5 2371.6   800.3 General Electric 1953
## 60  189.6 2759.9   888.9 General Electric 1954

Podemos utilizar a função a função plotmeans do pacote gplots para explorar graficamente os dados de Grunfeld:

gplots::plotmeans(invest ~ firm, main="Heterogeineidade entre empresas", 
                  n.label = FALSE, data = Grunfeld)

A inspeção do gráfico mostra que aparentemente há uma heterogeneidade substancial do investimento médio entre as empresas, entretanto, é necessário um teste formal para verificar o padrão observado.

gplots::plotmeans(invest ~ year, main="Heterogeineidade entre os anos", 
                  n.label = FALSE, data = Grunfeld)

O gráfico indica que aparentemente não há uma heterogeneidade substancial do investimento médio das empresas entre os anos até 1950, quando se verifica um aumento da variabilidade, mas, novamente, é necessário um teste formal para verificar o padrão observado.

Exemplo 2: Importanto e Estruturando Dados em Painel

O arquivo Data_AgricultureClimate.csv contém um painel envolvendo 568 municípios paulistas com informações anuais (1990 a 2014) para o valor da produção agrícola (em R$), área agrícola, temperatura média e precipitação total no ano.

Especificamente, o arquivo contém as seguintes variáveis:

• ano: 1990 a 2014
• munic: código do município no IBGE
• vtotal: valor total da produçao agropecuária (R$)
• areatotal: área agrícola no município
• tempano: temperatura média no ano
• precano: precipitacao total no ano

Estes dados foram utilizados no artigo de Maia at. al. (2018). Os autores
estimaram variações para o seguinte modelo para as variáveis em estudo:

\[ log(vtotal) = log(areatotal) + tempano + precano + erro \]

1. Importe o arquivo Data_AgricultureClimate.csv usando a função read_csv do pacote readr e salve o objeto com o nome agri. Visualize a estrutura dos dados com a função glimpse() do pacote dplyr e exiba estatísticas descritivas básicas do objeto agri usando a função interna summary()

painel <- "dados/Data_AgricultureClimate.csv"
agri <- read_csv(painel)
glimpse(agri)

## Rows: 14,200
## Columns: 8
## $ ano       <dbl> 1990, 1991, 1992, 1993, 1994, 1995, 1996, 1997, 1998, 1999, …
## $ vtotal    <dbl> 1264258.8, 737239.9, 1381728.0, 1547278.9, 1383899.7, 163200…
## $ areatotal <dbl> 1388970, 1235817, 1059597, 941193, 991584, 972081, 1063953, …
## $ tempano   <dbl> 22.86265, 23.21814, 22.73152, 23.19672, 23.60301, 23.55167, …
## $ tempanodp <dbl> 4.272846, 3.081071, 3.306779, 3.720999, 3.519357, 2.998808, …
## $ precano   <dbl> 1316.3407, 1414.5173, 1389.9879, 1319.8478, 1160.3283, 1492.…
## $ precanodp <dbl> 10.838446, 12.051548, 10.824051, 9.625086, 10.305651, 11.647…
## $ munic     <dbl> 3500105, 3500105, 3500105, 3500105, 3500105, 3500105, 350010…

summary(agri)

##       ano           vtotal           areatotal           tempano     
##  Min.   :1990   Min.   :       4   Min.   :       0   Min.   :14.34  
##  1st Qu.:1996   1st Qu.:  591117   1st Qu.:  227997   1st Qu.:20.53  
##  Median :2002   Median : 1598384   Median :  649044   Median :21.59  
##  Mean   :2002   Mean   : 2851582   Mean   : 1191138   Mean   :21.79  
##  3rd Qu.:2008   3rd Qu.: 3598698   3rd Qu.: 1487277   3rd Qu.:23.27  
##  Max.   :2014   Max.   :45904552   Max.   :14584086   Max.   :25.19  
##                 NA's   :1552       NA's   :375                       
##    tempanodp        precano         precanodp          munic        
##  Min.   :1.928   Min.   : 155.3   Min.   : 1.196   Min.   :3500105  
##  1st Qu.:2.954   1st Qu.:1226.1   1st Qu.: 8.290   1st Qu.:3514576  
##  Median :3.198   Median :1373.8   Median : 9.351   Median :3528552  
##  Mean   :3.193   Mean   :1406.4   Mean   : 9.531   Mean   :3528625  
##  3rd Qu.:3.432   3rd Qu.:1542.8   3rd Qu.:10.534   3rd Qu.:3542753  
##  Max.   :4.934   Max.   :4692.0   Max.   :26.780   Max.   :3557204  
## 

2. Declare o conjunto de dados como dados em painel usando a função pdata.frame() do pacote plm, nomeie o objeto criado como agri_pd e exiba a estrutura dos dados em painel usando a função pdim, também do pacote plm:

agri_pd <-  pdata.frame(agri, index = c("munic","ano"))
pdim(agri_pd)

## Balanced Panel: n = 568, T = 25, N = 14200

3. Exiba as primeiras 24 obeservações do objeto agri_pd

head(agri_pd, 48)

##               ano    vtotal areatotal  tempano tempanodp   precano precanodp
## 3500105-1990 1990 1264258.8   1388970 22.86265  4.272846 1316.3407 10.838446
## 3500105-1991 1991  737239.9   1235817 23.21814  3.081071 1414.5173 12.051548
## 3500105-1992 1992 1381728.0   1059597 22.73152  3.306779 1389.9879 10.824051
## 3500105-1993 1993 1547278.9    941193 23.19672  3.720999 1319.8478  9.625086
## 3500105-1994 1994 1383899.7    991584 23.60301  3.519357 1160.3283 10.305651
## 3500105-1995 1995 1632009.4    972081 23.55167  2.998808 1492.4040 11.647982
## 3500105-1996 1996 1910847.5   1063953 23.02875  3.409370 1152.8168  9.131358
## 3500105-1997 1997 1843234.7   1081377 23.31529  3.296005 1567.5603 10.897275
## 3500105-1998 1998 1812944.3   1046529 23.19375  3.455135 1389.1166  9.438434
## 3500105-1999 1999 2051107.5   1094841 23.47183  3.518524  980.5187  7.462520
## 3500105-2000 2000 1326583.7    969804 23.27866  3.855466 1268.0232  8.863539
## 3500105-2001 2001 2484606.3   1234827 23.63068  3.524630 1290.8452  9.092621
## 3500105-2002 2002 3081434.7   1445895 24.46877  3.400730 1026.9569  8.730663
## 3500105-2003 2003 2546356.2   1342341 23.49745  3.403117 1180.7192  8.541839
## 3500105-2004 2004 2503896.5   1463121 23.04267  3.649041 1224.7500  9.831753
## 3500105-2005 2005 2861892.0   1419561 23.41364  3.296626 1247.0122  9.437955
## 3500105-2006 2006 3008258.9   1273536 23.80656  3.295634 1511.0146 11.559646
## 3500105-2007 2007 3909081.1   1403127 24.15558  3.471190 1282.3023  9.389535
## 3500105-2008 2008 3452473.9   1560042 23.34564  3.164526 1083.7051  8.418116
## 3500105-2009 2009 3702484.3   1592217 23.35703  3.509578 1940.9464 14.676315
## 3500105-2010 2010 3326064.9   1571823 23.09775  3.375724 1312.7949 10.246465
## 3500105-2011 2011 3284135.7   1525095 22.95209  3.536957 1130.9613  8.458459
## 3500105-2012 2012 3603019.6   1514700 23.34242  3.645005 1269.6732 10.613371
## 3500105-2013 2013 3572662.4   1501038 22.53909  3.439130 1277.4219 11.539546
## 3500105-2014 2014 2756874.7   1415898 23.16971  3.339900 1117.6440  8.803157
## 3500204-1990 1990 1042861.2    607860 23.03726  3.703297 1314.6804  9.024720
## 3500204-1991 1991  970738.6    561330 22.99644  2.664130 1568.3458 11.120867
## 3500204-1992 1992 1671605.5    738936 23.07253  2.462134 1494.9482  9.391784
## 3500204-1993 1993 1784159.3    518859 23.35126  3.169225 1124.0356  7.474656
## 3500204-1994 1994 2412145.1    828729 23.92919  3.196905 1295.4416  9.031710
## 3500204-1995 1995 1035302.9    606870 23.72833  2.622399 1574.8358 10.609702
## 3500204-1996 1996  950784.9    552717 23.13837  3.316276 1099.6468  7.294775
## 3500204-1997 1997  729897.9    511038 23.66872  3.067165 1577.1214 10.494771
## 3500204-1998 1998 1052641.4    605088 23.70175  3.096906 1677.0027 11.680339
## 3500204-1999 1999 1157952.9    611919 23.61001  3.063488 1566.5552 12.579921
## 3500204-2000 2000  828805.4    498762 23.62397  3.344489 1397.3403  8.897663
## 3500204-2001 2001 3788505.3   1372536 23.82130  3.097323 1207.6119  8.130781
## 3500204-2002 2002 2011975.2    496188 24.65429  2.979756 1222.9100 11.267614
## 3500204-2003 2003 1955134.5    665280 23.65801  2.914628 1290.4623  8.640014
## 3500204-2004 2004 2685254.1    714681 23.39251  3.217841 1286.8070  9.379955
## 3500204-2005 2005 2385900.0    697752 23.98475  2.819513 1096.0276  8.250657
## 3500204-2006 2006 3211426.9    657360 23.70253  2.777337 1366.0578 10.178065
## 3500204-2007 2007 3089130.6    914265 24.26694  3.016799 1376.6647 10.531850
## 3500204-2008 2008 3530094.7    981387 23.37518  2.718268 1164.0335  8.075479
## 3500204-2009 2009 2865520.8    863775 23.04881  2.946778 1601.8338  9.154438
## 3500204-2010 2010 4999252.9    898128 23.44872  2.892584 1006.4249  8.112650
## 3500204-2011 2011 5397156.0    897138 22.82520  2.923734 1557.5957  9.926513
## 3500204-2012 2012 2640891.9    871299 23.54584  3.151972 1118.1704  7.209354
##                munic
## 3500105-1990 3500105
## 3500105-1991 3500105
## 3500105-1992 3500105
## 3500105-1993 3500105
## 3500105-1994 3500105
## 3500105-1995 3500105
## 3500105-1996 3500105
## 3500105-1997 3500105
## 3500105-1998 3500105
## 3500105-1999 3500105
## 3500105-2000 3500105
## 3500105-2001 3500105
## 3500105-2002 3500105
## 3500105-2003 3500105
## 3500105-2004 3500105
## 3500105-2005 3500105
## 3500105-2006 3500105
## 3500105-2007 3500105
## 3500105-2008 3500105
## 3500105-2009 3500105
## 3500105-2010 3500105
## 3500105-2011 3500105
## 3500105-2012 3500105
## 3500105-2013 3500105
## 3500105-2014 3500105
## 3500204-1990 3500204
## 3500204-1991 3500204
## 3500204-1992 3500204
## 3500204-1993 3500204
## 3500204-1994 3500204
## 3500204-1995 3500204
## 3500204-1996 3500204
## 3500204-1997 3500204
## 3500204-1998 3500204
## 3500204-1999 3500204
## 3500204-2000 3500204
## 3500204-2001 3500204
## 3500204-2002 3500204
## 3500204-2003 3500204
## 3500204-2004 3500204
## 3500204-2005 3500204
## 3500204-2006 3500204
## 3500204-2007 3500204
## 3500204-2008 3500204
## 3500204-2009 3500204
## 3500204-2010 3500204
## 3500204-2011 3500204
## 3500204-2012 3500204

gplots::plotmeans(vtotal ~ ano, main = "Heterogeineidade entre municipios", 
                  n.label = FALSE,
                  data = agri_pd)